/**
 * @author Ahammer
 * @date   2009 11
 * @update 2010 01 improved code, test image stack_0015.tif 5x enlarged to 930x1130 pixel: 10 times faster!
 * @update 2010 01 added progress bar support
 * @update 2010 03 added grey value images support DBC (Differential Box Counting), see Sarker Chauduri, IEEE Trans Syst Man Cyb 1994
 *                                                 RDBC Relative DBC, see Jin Ong Jayasooriah, Patt Rec Lett 16 1995   
 * @update 2011 01 added 16 bit support   
 * @update 2011 02 added regression data to table    
 * @update 2011 02 removing ModelListener extremely speeds up model configuration (adding columns, writing data,...)                                        
 */

/*
 * This file is part of Iqm.
 * Copyright (c) 2010-2011 Helmut Ahammer
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

package op;
import java.awt.RenderingHints;
import java.awt.image.DataBuffer;
import java.awt.image.Raster;
import java.awt.image.renderable.ParameterBlock;
import java.util.Vector;
import javax.media.jai.Histogram;
import javax.media.jai.ImageLayout;
import javax.media.jai.JAI;
import javax.media.jai.ParameterBlockJAI;
import javax.media.jai.PlanarImage;
import javax.media.jai.RenderedOp;
import javax.swing.JTable;
import javax.swing.table.DefaultTableModel;
import main.IqmTools;
import main.Manager;

/**
 * This is the main image processing class
 * There exist two approaches:
 * A user defined JAI operator is just called 
 * or:
 * The actual processing is implemented in this class 
 */
public class IqmFracBoxOperator{ 


  public IqmFracBoxOperator() {
	  //WARNING: Don't declare fields here
	  //Fields declared here aren't thread safe!   
  }
  /**
   * 
   * @param ParametrBlockJAI pb
   */
  public Object run(ParameterBlockJAI pbJAI){
	  	Object ob = null;
	  //Object ob = JAI.create("IqmFracBox".toLowerCase(), pb, null); //geht nicht
	  //daher hier direkt implementiert
	  	PlanarImage pi = (PlanarImage) pbJAI.getSource(0);
	  	int number     = pbJAI.getIntParameter("NumBoxes");
	  	int greyMode   = pbJAI.getIntParameter("GreyMode");
	  	int regStart   = pbJAI.getIntParameter("RegStart");
	  	int regEnd     = pbJAI.getIntParameter("RegEnd");
	  	boolean optShowPlot = false;
	  	boolean optDeleteExistingPlot = false;
	  	if (pbJAI.getIntParameter("ShowPlot") == 1 )          optShowPlot           = true;	  	
	  	if (pbJAI.getIntParameter("DeleteExistingPlot") == 1) optDeleteExistingPlot = true;	  	

	  	int numBands       = pi.getData().getNumBands();
	  	String type        = IqmTools.getImgTyp(pi);
	    double typeGreyMax = IqmTools.getImgTypeGreyMax(pi);
	  	boolean isBinary   = true;
//	  	int pixelSize = pi.getColorModel().getPixelSize();
//		String type = IqmTools.getCurrImgTyp(pixelSize, numBands);
	  	String imgName  = (String) pi.getProperty("image_name");
	  	String fileName = (String) pi.getProperty("file_name");
 
	  	//initialize table for Dimension data
		DefaultTableModel model = new DefaultTableModel();
		JTable jTable  = new JTable(model);	
		//adding a lot of columns would be very slow due to active model listener
		model.removeTableModelListener(jTable);
		model.addColumn("FileName");
		model.addColumn("ImageName"); 
		model.addColumn("Band");
		model.addColumn("RegStart"); 
		model.addColumn("RegEnd"); 
		for(int b = 0; b < pi.getNumBands(); b++){ //mehrere Bands
			model.addRow(new Object[]{fileName, imgName, b, regStart, regEnd});
		}
	  	double[][] totals  = new double[number][numBands];
	  	//double[] totalsMax = new double[numBands];	  	
		double[][] boxWidth = new double[number][numBands];
		double width = 1d;

		
	    // Set up the parameters for the Histogram object.
	    int[]    bins  = {(int)typeGreyMax+1, (int)typeGreyMax+1, (int)typeGreyMax+1}; // The number of bins  e.g. {256, 256, 256}
	    double[] lows  = {0.0D, 0.0D, 0.0D};                                           // The low incl.value  e.g. {0.0D, 0.0D, 0.0D}
	    double[] highs = {typeGreyMax+1, typeGreyMax+1, typeGreyMax+1};                // The high excl.value e.g. {256.0D, 256.0D, 256.0D}
	    // Create the parameter block.
	    ParameterBlock pbHisto = new ParameterBlock();
	    pbHisto.addSource(pi);   //Source image
	    pbHisto.add(null);       //ROI
	    pbHisto.add(1);          //Horizontal sampling
	    pbHisto.add(1);          //Vertical sampling
	    pbHisto.add(bins);       //Number of Bins
	    pbHisto.add(lows);       //Lowest inclusive values
	    pbHisto.add(highs);      //Highest exclusive values
	    // Perform the histogram operation.
		RenderedOp rOp = JAI.create("Histogram", pbHisto, null);
		// Retrieve the histogram data.
		Histogram histogram = (Histogram)rOp.getProperty("histogram");
		//System.out.println("IqmFracPyramidOperator: (int)typeGreyMax: " + (int)typeGreyMax);
	
		for (int b = 0; b < numBands; b++) {
			double totalGrey   = (double)histogram.getSubTotal(b, 1, (int)typeGreyMax);  //without 0!
			double totalBinary = (double)histogram.getSubTotal(b, (int)typeGreyMax, (int)typeGreyMax);
	    	if (totalBinary == totalGrey){
				isBinary = true;
			} else{
				isBinary = false;
			}
			if (isBinary) {
				totals[0][b] = totalBinary; //Binary Image
			}else{
				//First value
				//does not fit the regression line, too small!
				if (greyMode == 0){ //DBC
					totals[0][b] = 1.0*pi.getWidth()*pi.getHeight();  //Grey Image	(l-k+1) is always 1 because l=k	
				}
				if (greyMode == 1){ //RDBC
					//totals[0][b] = 0;  //Grey Image	 because u-b is always	zero
					totals[0][b] = 1.0*pi.getWidth()*pi.getHeight(); 
				}
			}
			
			//System.out.println("IqmFracBoxOperator: totals[0][b] "+ totals[0][b]);
			//totalsMax[b] = totals[0][b];
			//totals[0][b] = totals[0][b]/totalsMax[b]; 
			boxWidth[0][b] = 1d;
		}		 
		
	  	ParameterBlock pb = new ParameterBlock();
		pb.addSource(pi);
		pb.removeParameters();
		
	
		for (int n = 1; n < number; n++){
			int proz = (n+1)*100 / number;
   			Manager.setProgressBar(proz);
   			
//   			SwingUtilities.invokeLater( new Runnable() {
//				public void run() {
//					//GUI.parseProgressBar.setValue(i+1);
//					Manager.setProgressBar(proz);
//				}
//			});

			
			width = width * 2;
			pb.removeParameters();
			if (type == "RGB")    pb.add(DataBuffer.TYPE_BYTE);
			if (type == "8 bit")  pb.add(DataBuffer.TYPE_BYTE);
			if (type == "16 bit") pb.add(DataBuffer.TYPE_USHORT);
		
			ImageLayout layout = new ImageLayout();
			layout.setTileWidth((int) width);
			layout.setTileHeight((int) width);
			RenderingHints rh = new RenderingHints(JAI.KEY_IMAGE_LAYOUT, layout);

			pi = JAI.create("Format", pb, rh);

			Raster[] tileRaster = pi.getTiles();
			
			for (int b = 0; b < numBands; b++) totals[n][b] = 0d; 
			
			for (int r = 0; r < tileRaster.length; r++){ //tileRaster[]
				//slow version--------------------------------------------------------------------
				//only binary images!!!
				//and additionally problem with arbitrary image size (512 x 512 is ok, 800 x 800 is not ok)
//				WritableRaster wr = (WritableRaster) tileRaster[r];
//				WritableRaster newRaster = wr.createWritableTranslatedChild(0,0); 
//				BufferedImage bi = new BufferedImage(pi.getColorModel(), newRaster, false, null);
//				rOp = JAI.create("Histogram", bi);
//				histo = (Histogram)rOp.getProperty("histogram");
//				for (int b = 0; b < numBands; b++) {
//						double total = (double) histo.getSubTotal(b, 1, 255);
//						if (total > 0) totals[n][b] = totals[n][b] + 1d; 
//					////System.out.println("IqmFracBoxOperator: totals[n][b]  "+ totals[n][b] );
//					////totals[n][b] = totals[n][b]/totalsMax[b]; 
//					
//					////System.out.println("IqmFracBoxOperator: totals[n][b]  "+ totals[n][b] );
//				} //bands
			
				//fast version--------------------------------------------------------------------
				//test image stack_0015.tif 5x enlarged to 930x1130 pixel: 10 times faster!
				int minX = tileRaster[r].getMinX();
				int minY = tileRaster[r].getMinY();
				int tileWidth  = tileRaster[r].getWidth();
				int tileHeight = tileRaster[r].getHeight();
				for (int b = 0; b < numBands; b++) {
					boolean isGreaterZeroFound = false;
					int pixelMax = 0;                   //initial value
					int pixelMin = (int)typeGreyMax;    //initial value
					for(int x = 0; x < tileWidth; x++){
						for(int y = 0; y < tileHeight; y++){
							//System.out.println("IqmFracBoxOperator: b, x, y, " +b+"  " + x+ "  "+y);
							int pixel = tileRaster[r].getSample(minX + x, minY + y, b);
							if (isBinary){ //binary image
								if (pixel > 0) {
									totals[n][b] = totals[n][b] + 1d; 
									isGreaterZeroFound = true;
							}
							}else{ //grey image
								if (pixel > pixelMax) pixelMax = pixel;
								if (pixel < pixelMin) pixelMin = pixel;
							}
							if (isGreaterZeroFound) break;
							//System.out.println("IqmFracBoxOperator: totals[n][b]  "+ totals[n][b] );				
						} 
						if (isGreaterZeroFound) break;
					}
					if (!isBinary){
						double widthZ = width * typeGreyMax / ((pi.getWidth()+ pi.getHeight())/2.0);  //   epsilonZ/255 = epsilon/imageSize
						if (greyMode == 0){ //DBC
							double l = pixelMax / widthZ; //number of boxes in z direction
							double k = pixelMin / widthZ;
							l = Math.ceil(l);					
							k = Math.ceil(k);
							totals[n][b] = totals[n][b] + (l - k + 1);
						}	
						
						if (greyMode == 1){ //RDBC
							double diff = pixelMax - pixelMin;
							diff = Math.ceil(diff/widthZ);	
							totals[n][b] = totals[n][b] + diff;
						}
					}
				}
				
				//---------------------------------------------------------------------------------
								
			} // r tileRaster[]
			
			//normalization
			for (int b = 0; b < numBands; b++){	
				if (isBinary) {
					//do nothing
				}else{
					//totals[n][b] = totals[n][b]*width*width; //average volume
				}
			}
			
			for (int b = 0; b < numBands; b++){
				boxWidth[n][b] = (double) width;
			}

		} //n loop
		Manager.setProgressBar(0);

		//Log
		double[][] lnTotals  = new double[number][numBands];
		//double[][] lnNumbers = new double[number][numBands];
		double[][] lnBoxWidth = new double[number][numBands];
		for (int n = 0; n < number; n++){
			for(int b = 0; b < histogram.getNumBands(); b++){
				if (totals[n][b] == 0) totals[n][b] = Float.MIN_VALUE;  //damit logarithmus nicht undefiniert ist
				lnTotals[n][b]   = Math.log(totals[n][b]);
				lnBoxWidth[n][b] = Math.log(boxWidth[n][b]);
				//System.out.println("IqmFracBoxOperator: " );
			}
		}
		
		int numColumns = model.getColumnCount();
		model.addColumn("Db");
		model.addColumn("StdDev");
		model.addColumn("r2");
		model.addColumn("adjustet r2");
		for (int n = 0; n < number; n++) model.addColumn("DataX_"+(n+1));
		for (int n = 0; n < number; n++) model.addColumn("DataY_"+(n+1));
	
		boolean isLineVisible = false;
		for(int b = 0; b < histogram.getNumBands(); b++){ //several bands
			//Plot  //nur ein Band!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
			Vector<Double> dataX = new Vector<Double>();
			Vector<Double> dataY = new Vector<Double>();
			for (int n = 0; n < number; n++){	
					dataY.add(lnTotals[n][b]);
					dataX.add(lnBoxWidth[n][b]);
			}
			if (optShowPlot){
				IqmTools.displayRegressionPlotXY(dataX, dataY, isLineVisible,
												 imgName+"_Band"+b, "Box Dimension", "ln(Box Width)", "ln(Count)",
												 regStart, regEnd, optDeleteExistingPlot);
			}
			//double[] p = IqmTools.getLinearRegression(IqmTools.reverse(dataX), IqmTools.reverse(dataY), regStart, regEnd);
			double[] p = IqmTools.getLinearRegression(dataX, dataY, regStart, regEnd);

			double slope = p[1];
			double dim = 0.0;
			if (isBinary) {
				dim = -slope;
			}else{
				//dim = 3.0-(slope/2.0);
				dim = -slope;
			}
			model.setValueAt(dim,   b, numColumns); 
			model.setValueAt(p[3],  b, numColumns+1);	
			model.setValueAt(p[4],  b, numColumns+2);	
			model.setValueAt(p[5],  b, numColumns+3);
			//set regression data
			for (int n = 0; n < number; n++){	
				//model.setValueAt(lnBoxWidth[n][b],  b, numColumns+3+(n+1));
				model.setValueAt(boxWidth[n][b],  b, numColumns+3+(n+1));
			}
			for (int n = 0; n < number; n++){	
				//model.setValueAt(lnTotals[n][b],  b, numColumns+3+number+(n+1));
				model.setValueAt(totals[n][b],  b, numColumns+3+ number+(n+1));
			}
		}
		model.addTableModelListener(jTable);
		model.fireTableStructureChanged(); //this is mandatory because it updates the table
		jTable.getColumnModel().getColumn(2).setPreferredWidth(30);  //Band					  
		jTable.getColumnModel().getColumn(3).setPreferredWidth(30);  //RegStart				  
		jTable.getColumnModel().getColumn(4).setPreferredWidth(30);  //RegEnd					  	  
        return jTable;
	    //return ob;   //Optional wird ein Bild ausgegeben
  }
}
